Variational and diffusion Monte Carlo simulations of a hydrogen molecular ion in a spherical box
State Key Laboratory of Superhard Materials, College of Physics, Jilin University, Changchun 130012, China
† Corresponding author. E-mail:
baokuo@jlu.edu.cn cuitian@jlu.edu.cn
1. IntroductionIn recent years, studies focused on the change in atomic and molecular system properties during confinement inside penetrable or impenetrable boundaries. Atoms or molecules are known to exhibit different electrical and structural behaviors when squeezed into a small space. Therefore, determining how these changes occur is important. With the emergence of modern nanostructured materials, such as carbon nanotubes, which are ideal containers for molecular insertion and storage, research on confined molecules has become increasingly meaningful. Researchers have recently succeeded in inserting atoms and molecules into fullerene cages and derivatives in experiments.[1] Additionally, confined molecules can be applied to the study of high-pressure materials with molecular hydrogen.[2] The exploratory model must be supplemented to better understand the basic change mechanism in the electronic and structural properties of the restricted molecules and atoms.
The hydrogen molecular ion
, with only one electron and two protons, is one of the simplest molecular ion systems in the universe. In recent years, many researchers have studied the properties of the hydrogen molecular ion and its isotopes, such as the cold trapped molecular ion[3,4] and the effect of
in the interaction of H2 with transition metals.[5] Recently, Silva et al.[6] calculated the energy and polarizability of
confined to an ellipsoidal box in the ground and the first excited states with the variational method. The molecules are confined to an ellipsoidal box to study the properties of the molecules, which can also be called molecule-in-a-box model. This model is simple but effective. The molecule-in-a-box model was first proposed[7] to calculate the kinetic energy of one hydrogen atom. Over the years, this model has been extensively used, especially in small molecular systems, such as hydrogen and helium systems. The ground-state energy of the hydrogen molecule and molecular ion in the ellipsoidal box was reported to fix the protons in the focus position of the box by LeSar and Herschbach[8,9] using a five-term James–Coolidge variational function and by Pang[10] using quantum Monte Carlo simulation. Although previous researchers acquired some meaningful results,[11–13] they ignored the zero-point motion of protons. Cruz et al. first attempted the creation of protons away from the special point under the hydrogen molecular ion system[14] and then applied this change to the system.[15,16] They relaxed the position of protons to allow free movement on the semi-major axis of the ellipsoid, resulting the change in the equilibrium distance between protons and the total ground-state energy of the hydrogen molecules corresponding to fixed size and shape. Assuming that the protons can move freely in the ellipsoid instead of only moving on a semi-major axis, we continue to make improvements to this model in this work.
Given that the model is a typical quantum atom system, theoretical study on hydrogen is rather difficult. The researchers managed to develop three main theoretical methods, namely, density functional theory,[17–19] molecular dynamics simulation,[20,21] and Monte Carlo (MC) simulation.[22–24] For simplicity of treating the dissociation of hydrogen[25] and the zero-point motion of protons, we choose MC simulation to study the system in this work. Considering the degrees of freedom of electrons and protons, we extend the previous molecule-in-a-box model, and study the hydrogen molecular ion in a spherical box as a three-body problem. The variational Monte Carlo (VMC) and diffusion Monte Carlo (DMC) methods are used to accurately solve the three-body problem under the confinement condition, and their calculation precisions are enough to reach the required accuracy.
This paper is organized in the following manner. Section 2 provides a detailed theoretical description, including the theoretical model, trial wave function, and Markov Chain. Section 3 shows our results on the properties of the ground state and the comparison of these results with previous studies[6,8] to analyze differences in findings. The last section concludes the paper and presents further discussion.
2. Theory2.1. The theoretical modelThe original idea of the molecule-in-a-box model is that one can place
into an ellipsoidal box, where the protons are fixed at the focal points of the semi-major axis, and the electron can move freely in the ellipsoid. Given that the protons are fixed, this model simulates the properties under different constraint conditions by changing the volume or shape of the ellipsoidal box, such as the ground-state energy, equilibrium bond length, and kinetic energy of electron. The model has been applied to simulate the pressure effect of atoms, molecules, and other bound electron systems.[6,8,10,26–28] Although the molecule-in-a-box model works well when dealing with various properties of matter, it overestimates the effect of compression in most cases. We use a revised version of the molecule-in-a-box model in the present work. By using this configuration as the starting point, the electron and proton can move freely inside the box to follow the Markov chain.
For the molecule-in-a-box model, the boundary is defined by
where
a and
b are the semi-major and semi-minor axes, respectively.
The Hamiltonian of this system can be written as
where
Te (
Tp) represents the kinetic energy of electron (protons) and
V(
R) is the Coulomb potential energy of this system, which can be written as
In this paper, subscript 1 represents the only electron, and subscripts 2 and 3 represent the protons. The previous molecule-in-a-box model calculation did not consider the motion of protons and did not have the term
Tp. Given that only two protons and one electron exist in the hydrogen molecular ion, the potential energy only contains two items, namely, the potential energy between the electron and protons and the potential energy between the two protons. The setting of atomic units used in this work is
.
The Schrödinger equation can be written as follows:
where
is the chosen trial wave function, and
E is the ground-state energy of
. In this work, we use VMC and DMC methods to simulate
in a spheroidal box as a three-body problem.
2.2. The trial wave functionThe quality of the trial wave function can affect the calculation speed and results in VMC and DMC simulations. Therefore, selecting the trial wave function is crucial. The chosen trial wave function in this work contains the molecular orbital wave function of electron, the movement of protons, and the boundary conditions. This function is written in the following form:
where
is the molecular orbital wave function, which is written as the linear combination of atomic orbital
[10]
denotes the movement of protons
[29,30] and is written as
and
describes the boundary conditions caused by the spheroidal box
[10] as follows:
The wave function satisfies the zero value when the electron or proton is outside of the spheroidal box. In the above mentioned trial wave function,
α,
C1, and
C2 are the variational parameters. In the process of VMC simulation, we can first optimize the chosen trial wave function and ground-state energy of
to determine the variational parameters. Simultaneously, the semi-major axis and corresponding semi-minor axis are adjusted to determine their values at the lowest energy. On the basis of the VMC simulation, we use the DMC method to reconfirm the values of semi-major and semi-minor axes for a given volume and accurately calculate the related physical quantities of this system.
2.3. Markov chainIn the calculation process, we explore the movement of the proton and electron by using the Markov chain in the MC simulation. The movement of the proton and electron can be respectively written in the following forms:
where
Xe (
) and
Xp (
) are the vectors in terms of electron and proton positions, respectively;
δ is a random number between 0 and 1; and
and
represent the maximum distances of electron and proton movement each time, respectively. The preceding formula ensures that the electron (proton) randomly moves between
(
) and 0.5
(
) in the ellipsoid.
Xe and
Xp are the positions of electron and proton at this moment, respectively, which are subjected to a dynamic rule process (Eqs. (
11) and (
12)) to randomly generate the next positions
and
. This phenomenon is the concept of the Markov chain, which is one of the important factors of stochastic simulation.
3. ResultsWe calculate the ground-state energy of
in different volumes using two MC simulation methods. Figure 1 shows the relationship between the ground-state energy and the box volume, where the energy increases with the decrease in volume. The energy slightly changes with the decrease in volume when the volume is larger than 200
(a0 is the Bohr radius). However, the energy rapidly changes when the volume is smaller than 200
mainly due to the increase in pressure and system energy as the volume decreases. Under the same conditions, our results are similar to those of LeSar and Silva but slightly larger, mainly because we consider the motion of the proton and include its kinetic energy.
Through the change in the ground-state energy of this system along with the bound volume, we can obtain the pressure corresponding to the given volume in accordance with the equation
. Figure 2 provides the equation of states (EOS) of the
system. In Fig. 2, the calculated pressures are larger than those of LeSar for the same density. When the density is small, the calculated pressures and those of LeSar are almost the same. However, with the increase in density, our results are much higher than their calculated results. Given the limitation of this model, the calculations of this study and those of LeSar are slightly larger than the real values. The largest difference between our calculation and that of others is that the protons are not fixed in the special points, such as the focuses of the spheroidal box. Thus, the zero-point motion of the protons becomes evident under large density.
The major difference between the present work and the previous studies is the assumption that the protons and electron can move freely in the same way in the simulation process as emphasized earlier. Therefore, the calculations of the distances between protons under different pressures are crucial. Figure 3 presents the calculation results for the mobile model and those obtained by the previous fixed model. In Fig. 3, the distances calculated in the mobile model are larger than those calculated in the fixed model, but the trends of the two curves are the same. The distances between the protons decrease as the pressures increase. On our calculated curve, the distance is observed to have a slight increase at a pressure of approximately 0.4 GPa. This phenomenon may be related to the reduction in the electron density in the box.[16] Notably, the calculated distance in the volume of 879
is experimentally similar to the bond length of the hydrogen molecular ion. Our result is 2.09a0, and the experimental result is 2.08a0. This finding also proves that the chosen trial wave function is reasonable and our calculations are accurate.
We calculate the average
in the VMC simulation to better understand the motion of electron and protons and the effects of considering the motion of protons. Figure 4 shows that the average
in the mobile model is larger than that in the fixed model, and the trends of the two curves are the same. The average
increases with the increase of the semi-major axis. However,
in the mobile model increases rapidly with the semi-major axis of less than 4.0 a.u., and the growth of
is slowly or even basically unchanged when the semi-major axis exceeds 4.0 a.u. In the fixed model, the turning point is around 3.0 a.u., and
still increases slowly when the semi-major axis is larger than 3.0 a.u., indicating that the motion of protons can affect the motion of the electron.Considering the motion of the proton, the change trend of
in the mobile model is obviously different from that in the fixed model, which is similar to a certain phase transition. However,
is just a calculating parameter, which is not directly related with measurements. The discontinuous
does not necessarily indicate a traditional phase transition, which generally means that due to the changes in temperature, pressure, or other external conditions, certain properties of the medium usually change discontinuously. In order to observe the impact of proton motion on electron motion and verify if there is a phase transition, we calculate the corresponding kinetic energy of the electron (Te) at different volumes in the two models and the relationship between Te and 1/V. Figure 5(a) demonstrates that the kinetic energy decreases with the increase in box volume, and Te in the mobile model is larger than that in the fixed model. We calculate Te under the two models in detail (see Table 1) and find that the kinetic energy increases by about 5%. Combining Figs. 4 and 5(a) for the same model, the larger the average
is, the smaller the Te is. In Fig. 5(b), the variation trends of the curves in the two models are basically the same, and Te and 1/V are basically proportional to each other. We think that there is no obvious phase transition. 6
Finally, we discuss the relationship of the ratio of the kinetic energy of proton and electron (Tp/Te) and pressure. Overall, Tp/Te exponentially decreases with increasing pressure. Tp/Te is small and always of the order of a few thousandths in magnitude, which is of the same order of magnitude as the ratio of mass of electron and proton. Given that the proton kinetic energy accounts for a small proportion of the total kinetic energy, this energy may be ignored if the calculation accuracy is extremely low, which is related to the specific model and the discussed problem. This result provides a valuable reference for the calculation method, which fails to fully estimate the influence of proton motion on the possible error.
4. ConclusionIn summary, we employ the improved molecule-in-a-box model which is simple but effective in calculating the ground-state properties of the hydrogen molecular ion under high pressure. We consider the zero-point motion of protons and regard the hydrogen molecular ion as a three-body system. Our results are larger than the previous results in the fixed model due to the influence of proton motion. In this work, the results show that the total energy of the system decreases with the increase in volume and the distance between protons decreases with the increase in pressure. We also calculate the EOS of the hydrogen molecular ion, and our results are similar to those of our predecessors under low density (low pressure), but higher than theirs in the high density range. We quantitatively study the effect of the zero-point motion of protons on electrons. The kinetic energy of electrons in the mobile model increases by 5% compared with that in the fixed model. Finally, the study shows that the kinetic energy of the protons is small compared with that of the electron, which is only a few thousandths of Te. Therefore, Tp can be neglected when the calculation accuracy is low.
